// Based on:
//
// Jupiter.java -- calculate positions of Jupiter's moons.
//
// Equations come from Jean Meeus, Astronomical Formulae for Calculators.
//
//	Copyright 1996 by Akkana Peck.
// You may use or distribute this code under the terms of the GPL.
//
// $Rev: 54 $
// $LastChangedDate: 2009-07-15 14:22:06 +0000 (Wed, 15 Jul 2009) $

package Main;

/**
 *
 * @author moreto
 */

public class Jupiter
{
    final static int NUM_MOONS = 4;

    double d;		// days since 1899 Dec 31 12h ET

    // Angles of each of the Galilean satellites, in radians,
    // expressed relative to each satellite's inferior conjunction:
    double moonAngles[];
    // And their distances from the planet:
    double moonDist[];
    // And their abbreviations:
    String moonLabels[];

    double delta;	// Earth-Jupiter distance
    double De;		// planetocentric ang. dist. of earth from jup. equator
    double G, H;

    // latitudes of systems I and II:
    double lambda1;
    double lambda2;

    // variables we may want later for moon calcs:
    double psi;

    boolean topFlag = false;

    public Jupiter(StarDate _date, boolean flag){
        moonAngles = new double[NUM_MOONS];
        moonDist = new double[NUM_MOONS];
        moonLabels = new String[NUM_MOONS];
        moonLabels[0] = "Io";
        moonLabels[1] = "Eu";
        moonLabels[2] = "Ga";
        moonLabels[3] = "Ca";

        topFlag = flag;
        setDate(_date);
    }

    int getNumMoons(){
        return NUM_MOONS;
    }

    void setTopFlag(boolean flag){
        topFlag = flag;
    }
    //
    // Convert an angle (in radians) so that it's between 0 and 2*PI:
    //
    double angle(double a){
        if (a < 10000)
            return oangle(a);
        a = a - 2.*Math.PI * (int)(a / 2. / Math.PI);
        if (a < 0)
    	    a += 2*Math.PI;
        return a;
    }

    double oangle(double a)
    {
        while (a > 2 * Math.PI)
            a -= 2. * Math.PI;
        while (a < 0)
            a += 2. * Math.PI;
        return a;
    }

    void setDate(StarDate jd)
    {
        // Calculate the position of Jupiter's central meridian:

        // First, get the number of days since 1899 Dec 31 12h ET.
        d = jd.getJulianDate() - 2415020.5;     // days since 1899 Dec 31 12h ET

        // Argument for the long-period term in the motion of Jupiter:
        double V = angle((134.63 + .00111587 * d) * Math.PI / 180);

        // Mean anomalies of Earth and Jupiter:
        double M = angle((358.476 + .9856003 * d) * Math.PI / 180);
        double N = angle((225.328 + .0830853 * d + .33 * Math.sin(V))
			 * Math.PI / 180);

        // Diff between the mean heliocentric longitudes of Earth & Jupiter:
        double J = angle((221.647 + .9025179 * d - .33 * Math.sin(V))
			 * Math.PI / 180);

        // Equations of the center of Earth and Jupiter:
        double A = angle((1.916 * Math.sin(M) + .020 * Math.sin(2*M))
        	 	 * Math.PI / 180);
        double B = angle((5.552 * Math.sin(N) + .167 * Math.sin(2*N))
        		 * Math.PI / 180);

        double K = angle(J + A - B);

        // Distances are specified in AU:
        // Radius vector of the earth:
        double R = 1.00014 - .01672 * Math.cos(M) - .00014 * Math.cos(2*M);
        // Radius vector of Jupiter:
        double r = 5.20867 - .25192 * Math.cos(N) - .00610 * Math.cos(2*N);

        // Earth-Jupiter distance:
        delta = Math.sqrt(r*r + R*R - 2*r*R*Math.cos(K));

        // Phase angle of Jupiter (always btw. -12 and 12 degrees):
        psi = mMath.asin(R / delta * Math.sin(K));

        // Longitude of system 1:
        lambda1 = angle((268.28 * 877.8169088 * (d - delta / 173))
			* Math.PI / 180 + psi - B);
        // Longitude of system 2:
        lambda2 = angle((290.28 + 870.1869088 * (d - delta / 173))
			 * Math.PI / 180 + psi - B);

        // calculate the angles of each of the satellites:
        moonAngles[0] = angle((84.5506 + 203.4058630 * (d - delta / 173))
			      * Math.PI / 180
			      + psi - B);
        moonAngles[1] = angle((41.5015 + 101.2916323 * (d - delta / 173))
			      * Math.PI / 180
			      + psi - B);
        moonAngles[2] = angle((109.9770 + 50.2345169 * (d - delta / 173))
			      * Math.PI / 180
			      + psi - B);
        moonAngles[3] = oangle((176.3586 + 21.4879802 * (d - delta / 173))
			      * Math.PI / 180
			      + psi - B);

        // and the planetocentric angular distance of the earth
        // from the equator of Jupiter:
        double lambda = angle((238.05 + .083091 * d + .33 * Math.sin(V))
                              * Math.PI / 180 + B);
        De = ((3.07 * Math.sin(lambda + 44.5 * Math.PI / 180)
	       - 2.15 * Math.sin(psi) * Math.cos(lambda - 24.*Math.PI/180)
	       - 1.31 * (r - delta) / delta
	       * Math.sin(lambda - 99.4 * Math.PI / 180))
	      * Math.PI / 180);

        if (topFlag){
            De = 90;
        }

        G = angle((187.3 + 50.310674 * (d - delta / 173)) * Math.PI / 180);
        H = angle((311.1 + 21.569229 * (d - delta / 173)) * Math.PI / 180);

        // Calculate the distances before any corrections are applied:
        moonDist[0] = 5.9061 -
		       .0244 * Math.cos(2 * (moonAngles[0] - moonAngles[1]));
        moonDist[1] = 9.3972 -
		       .0889 * Math.cos(2 * (moonAngles[1] - moonAngles[2]));
        moonDist[2] = 14.9894 - .0227 * Math.cos(G);
        moonDist[3] = 26.3649 - .1944 * Math.cos(H);

        // apply some first-order correction terms to the angles:
    	moonAngles[0] = angle(moonAngles[0] +
			      Math.sin(2 * (moonAngles[0] - moonAngles[1]))
			      * .472 * Math.PI / 180);
    	moonAngles[1] = angle(moonAngles[1] +
			      Math.sin(2 * (moonAngles[1] - moonAngles[2]))
			      * 1.073 * Math.PI / 180);
       	moonAngles[2] = angle(moonAngles[2] +
			      Math.sin(G) * .174 * Math.PI / 180);
    	moonAngles[3] = angle(moonAngles[3] +
			      Math.sin(H) * .845 * Math.PI / 180);
    }

    //
    // Routines to get distances, angles, or X and Y coordinates
    // for each of the moons.
    // These assume that moon positions have already been calculated!
    //
    double getMoonDist(int whichmoon){
        return moonDist[whichmoon];
    }
    double getMoonAngle(int whichmoon){
        return moonAngles[whichmoon];
    }
    String getMoonLabel(int whichmoon){
        return moonLabels[whichmoon];
    }

    XYCoord getMoonXY(int whichmoon)
    {
	double r = getMoonDist(whichmoon);

    	XYCoord coord = new XYCoord();

	// See whether the moon is hidden behind the planet:
	coord.x = r * Math.sin(moonAngles[whichmoon]);
	coord.y = -r * Math.cos(moonAngles[whichmoon]) * Math.sin(De);

	if (coord.x < 1. && coord.x > -1.
	    && moonAngles[whichmoon] > Math.PI * .5
	    && moonAngles[whichmoon] < Math.PI * 1.5)
	{
	    coord.x = coord.y = Double.POSITIVE_INFINITY;
	    // NOTE!  Can't use NaN, because comparisons with !=
	    // don't work against NaN.
	}
    	return coord;
    }

    //
    // Moon shadows
    //
    XYCoord getMoonShadowXY(int whichmoon)
    {
	double r = getMoonDist(whichmoon);
	double moonSunAngle = moonAngles[whichmoon] - psi;

    	XYCoord coord = new XYCoord();
        coord.x = r * Math.sin(moonSunAngle);
        // This Y coord isn't right ... need to derive the right eqn:
        coord.y = -r * Math.cos(moonSunAngle) * Math.sin(De);

	// See whether the shadow is on the far side of the planet:
	if (coord.x > -1. && coord.x < 1.
	    && moonSunAngle > Math.PI * .5 && moonSunAngle < Math.PI * 1.5)
	{
	    coord.x = coord.y = Double.POSITIVE_INFINITY;
	}
        else if (coord.x < -1. || coord.x > 1.)
        {
            // shadow isn't hitting the planet right now
            // Some day, ought to check for moons eclipsing other moons
	    coord.x = coord.y = Double.POSITIVE_INFINITY;
        }

    	return coord;
    }

    //
    // The Great Red Spot, currently at longitude 61 in system II
    //
    XYCoord getRedSpotXY(int spot_in_deg)
    {
	double spotlong = angle(lambda2 - spot_in_deg*Math.PI/180);

    	XYCoord coord = new XYCoord();

	// See if the spot is visible:
	if (spotlong > Math.PI * .5 && spotlong < Math.PI * 1.5) {
	    coord.x = coord.y = Double.POSITIVE_INFINITY;
	} else {
	    coord.x = Math.sin(spotlong);
	    coord.y = .42;	// completely random wild-assed guess
	}

    	return coord;
    }

    //
    // You might also want to get the location of some arbitrary
    // other position on the planet, e.g. the Great Northern Spot.
    //
    double getJovianPointX(double long_in_deg, int systm){
        double lambda = (systm == 1 ? lambda1 : lambda2);
        double longInRad = angle(lambda - long_in_deg*Math.PI/180);

        // See if the point is visible:
        if (longInRad > Math.PI * .5 && longInRad < Math.PI * 1.5) {
            return Double.POSITIVE_INFINITY;
        } else {
            return Math.sin(longInRad);
        }
    }
}
